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Impurities cause radiation losses and plasma dilution, and in stellarator plas¬ 
mas the neoclassical ambipolar radial electric field is often unfavorable for avoiding 
strong impurity peaking. In this work we use a new continuum drift-kinetic solver, 
the SFINCS code (the Stellarator Fokker-Planck Iterative Neoclassical Conserva¬ 
tive Solver) [M. Landreman et al., Phys. Plasmas 21 (2014) 042503] which employs 
the full linearized Fokker-Planck-Landau operator, to calculate neoclassical impu¬ 
rity transport coefficients for a Wendelstein 7-X (W7-X) magnetic configuration. 

We compare SFINCS calculations with theoretical asymptotes in the high collision- 
ality limit. We observe and explain a l/i'-scaling of the inter-species radial transport 
coefficient at low collisionality, arising due to the field term in the inter-species col¬ 
lision operator, and which is not found with simplified collision models even when 
momentum correction is applied. However, this type of scaling disappears if a radial 
electric field is present. We also use SFINCS to analyze how the impurity con¬ 
tent affects the neoclassical impurity dynamics and the bootstrap current. We show 
that a change in plasma effective charge Z e Q of order unity can affect the bootstrap 
current enough to cause a deviation in the divertor strike point locations. 
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I. INTRODUCTION 


3D plasma confinement concepts have an advantage over the tokamak, as they offer 
the potential of steady state plasma operation with no need for current drive [1, 2]. For 
steady state operation, impurity accumulation has to be avoided, since impurities cause 
plasma dilution, radiation losses and can lead to pulse termination by radiation collapse. 
The avoidance of impurity accumulation under relevant conditions is one of the most 
crucial tests for the capability of steady state 3D systems. 

In stcllarators, particles can be trapped in helical magnetic wells and escape the 
plasma even in the absence of collisions. Thus the neoclassical transport is typically 
considerably larger than in axisymmetric configurations. In fact, at low collisionality 
(such as in a hot plasma core) neoclassical transport could be expected to dominate 
over the turbulent transport because of the l/i/-transport behavior, v being the collision 
frequency [2]. The flux-surface-averaged radial neoclassical particle flux of species a can 
be written as 


(r o • V^) = Ij d 3 vfaiv da • wA = -n a ^2 
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where n& is the density, Tb is the temperature and e& = Zbe is the charge of species b with e 
being the proton charge, and the sum is taken over all plasma species. The brackets (...) 
denote a flux-surface average. The Dj b :s are coefficients of the transport matrix, \ da is the 
cross-held drift velocity and f a i = f a — fMa is the departure from the Maxwellian part of 
the distribution function of species a. ij) is a flux function representing a radial coordinate 
(often chosen to be the toroidal magnetic flux) and $ is the electrostatic potential which 
relates to the radial electric held by E r = —chh/dr where r is the effective radius. The 
inter-species terms ( b ^ a) are due to friction along the magnetic held between the different 
species [2]. Moreover, in all stcllarator collisionality regimes D% a is positive which implies 
that the temperature gradient drives an outward particle hux, tending to generate hollow 
density prohles. 

A consequence of collisionless trajectories not necessarily being confined is that 
different plasma species can have different radial transport rates. This results in a radial 
electric held E r to restore ambipolarity, which can be determined without knowledge of 
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the turbulent transport since the radial neoclassical current is 1/p* = L/pi (p* being 
the ion gyro radius and L a typical macroscopic scale length) larger than the radial 
turbulent current unless E r is just right [3]. We note that the transport coefficients, 
Dj b , in Eq. (1) depend on the value of E r . The ambipolarity condition of the particle 
fluxes determining E r can have multiple roots, depending on the plasma parameters and 
the magnetic configuration. In the standard situation for stellarators, the neoclassical 
ambipolar radial electric field points radially inwards (i.e. a negative E r ), referred to as 
the ion root regime (in the electron root regime E r is instead positive). This electric held 
tends to cause impurity accumulation, which is particularly strong for heavy impurities 
whose charge numbers Z arc large. The electrostatic drive for impurity accumulation is not 
present in tokamaks or quasi-symmetric stellarators, where intrinsic ambipolarity implies 
that radial transport rates of each species must be independent of the radial electric held 
to high precision, and accumulation is merely caused by impurity-ion friction. 

In a single-impurity species plasma, we denote the species by the subscripts e 
(electrons), i (bulk hydrogen ions) and z (impurities) respectively. For a trace impu¬ 
rity approximation (when Zn z /n e <C 1) in a standard ion root plasma, the ambipolar 
electric held is essentially determined from the condition (Tj (E r ) ■ \7p)) = 0 because 
(Tj (E r = 0) • V-^) (T e • Vr/) [4]. Assuming T z = Ti = T and neglecting the inter¬ 
species coefficients in Eq. (1), it is possible to express the radial neoclassical impurity hux 
as 
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where the term containing the factor Z appears from substituting the ambipolarity con¬ 
dition for the bulk ions. Df z represents the diffusive part connected to gradients in n z , 
while V z are the convective terms related to gradients of the bulk species. We note that 
Eq. (2) is also valid for an axisymmetric device, but the way to derive it differs from the 
way to do it for a stellarator because the particle hux in Eq. (1) is independent of the 
radial electric held in axisymmetry. For tokamaks the coefficient Dip/D\ l in Eq. (2) can 
be negative indicating impurity screening. In stellarators however, earlier theory predicts 
it to be always positive and support impurity accumulation when the temperature prohle 
is peaked [4], Typically a transient increase in impurity concentration is found, due to 
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the inward convective impurity transport which drives impurity accumulation, until it 
balances the outward diffusive transport and (r 2 • V^>) = 0. This results in a peaking 
of the impurity density profile with respect to the main ion profile, depending on the 
relative strength of the impurity pinch V z /Dl z . This stationarity in the impurity profile 
is expected ultimately, even if the impurity core confinement time is large. 

Although neoclassical predictions often point towards impurity accumulation, there 
are experimental scenarios in stcllarators where impurity accumulation has been avoided, 
e.g. in certain low-density scenarios at W7-AS and LHD with an outward radial electric 
held, or through the application of purification mechanisms such as radiofrequency heat¬ 
ing. In LHD an extremely hollow profile of carbon impurity has been observed, referred to 
as an “impurity hole” [5, 6]. In W7-AS also a “high density H-mode” with low impurity 
confinement times has been discovered, accessible only through neutral beam injection 
[7]. It is also possible that turbulent transport could significantly mitigate the neoclassi¬ 
cal impurity accumulation. To enable the stcllarator concept as a fusion reactor candidate 
with high pressure plasmas, the search for favorable scenarios capable of avoiding strong 
impurity accumulation is important. 

In neoclassical theory, transport processes are usually assumed to be radially lo¬ 
cal and described by the linearized drift kinetic equation for the first order distribution 
function f a \ in p*. To date, stcllarator neoclassical calculations have predominantly been 
performed using simplified models for collisions. The most accurate linear operator avail¬ 
able is the linearized Fokker-Planck-Landau operator [8, 9] (also given in Appendix A), 
which is used in a variety of axisymmetric calculations. However because of the extra 
dimension in stcllarators, due to the lack of toroidal symmetry, stellarator calculations 
are more challenging and often only pitch-angle scattering collisions are retained (e.g. see 
Ref. [10]). This implies that coupling in the energy dimension is eliminated, and that 
momentum is generally not conserved. Momentum correction methods exist [11-13] for 
post-processing pitch-angle scattering results, but these methods are not equivalent to 
using the full linearized Fokker-Planck-Landau operator. The difference between calcu¬ 
lations with momentum-corrected pitch-angle scattering and full Fokker-Planck-Landau 
collisions could be expected to be particularly important for ion-impurity and impurity-ion 
collisions since the mass ratio is neither very large nor very small. 
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In this work we study neoclassical impurity transport in stcllarators using the con¬ 
tinuum code SFINCS (the Stcllarator Fokker-Planck Iterative Neoclassical Conservative 
Solver), described in [14]. The code solves the radially local 4D drift-kinetic equation, 
retaining coupling in four of the independent phase space variables (two spatial and two in 
velocity). The code permits an arbitrary number of species, and it includes the linearized 
Fokker-Planck-Landau operator for self- and inter-species collisions, with no expansion 
made in mass ratio. The present numerical implementation of this operator in the code is 
detailed in Appendix A. Our study will be restricted to a hydrogen plasma with one single 
impurity species. We emphasize that we will focus solely on neoclassical transport in this 
work, and stellarator turbulent transport is an area where much is still to be explored [2], 
Moreover, Ref. [15] shows that the neoclassical impurity transport can be strongly affected 
by the variation of the electrostatic potential on a flux-surface, <!>! = <f> — (<f>), which can 
be large enough to affect impurity species of high charge. Although SFINCS calculates 
4>i, this effect is not included in the calculations presented here. 

The remainder of the paper is organized as follows. In Sec. II we use SFINCS to 
calculate neoclassical transport coefficients for the impurity particle flux, in a single- 
impurity-species hydrogen W7-X plasma. We discuss the importance of the full linearized 
Fokker-Planck-Landau operator, by comparing to results from pitch-angle scattering cal¬ 
culations where momentum correction is applied afterwards. Furthermore, we investigate 
in which collisionality regimes impurity screening from the pressure gradients can be 
expected. How the impurity content affects the neoclassical impurity dynamics and the 
bootstrap current in a non-axisymmetric device is then explored with SFINCS in Sec. III. 
In Sec. IV we summarize the results and conclude. 
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II. IMPURITY TRANSPORT COEFFICIENTS 


In this section we will calculate neoclassical transport coefficients L zb k for the impu¬ 
rity particle flux, which is written as a linear combination of the thermodynamic forces 


1 dn z Ze d<3> 3 dT z 
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( 3 ) 


where 2nip from here on is specifically the toroidal magnetic flux. The motivation for our 
choice of the thermodynamic forces in Eq. (3) stems from Ref. [14], where the transport 
matrix for a single species drift-kinetic system of equations becomes Onsager symmetric 
if E r — 0 and the forces are defined in this form. We assume a hydrogen plasma with 
a single impurity species present where the impurities and the main ions are in thermal 
equilibrium, T z = Tj = T. We consider a C 6+ impurity, as this species is expected to be 
the dominant impurity in W7-X. Note that we neglect the impurity-electron collisions, 
since their effect on the collisional impurity transport is negligible whenever Z 2 n z /ni ^> 
y 7 m e /mi , which is practically always the case in reality [16]. 

It is possible to express the impurity and ion fluxes in terms of the transport matrix 
L°j b defined as follows: 
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The normalization is similar to Eq. (40) of [14], and implies that the matrix elements are 
dimensionless. Moreover, the matrix elements in Eq. (4) depend on E r , and it is possible 
to show that L°j b k (E r ) = L b k - (—E r ) when = T z . For a stellarator-symmetric magnetic 
geometry the elements are independent of the sign of E r , L? k (. E r ) = L^ k (—E r ). Hence, 
written in this form the matrix exhibits Onsager symmetry, i.e. L°j k = L b k -. We note 
that the Onsager symmetry is true for both the “full particle trajectories” described by 
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Eq. (17) and the “DKES particle trajectories” described by Eq. (18) in [14]. If we were 
to include additional matrix elements in Eq. (4) corresponding to bootstrap current and 
Ware pinch, the Onsager symmetry would generally not be fulhlled for the “full particle 
trajectories” when E r ^ 0. In this work the “full particle trajectories” is the default in the 
SFINCS calculations. However, we will also present results from SFINCS calculations 
with “DKES particle trajectories” for E r ^ 0 (when E r = 0 the two different models for 
the particle trajectories are equal). 

In this work we focus on the impurity particle transport, and write 

< r '' Vl » = zSk ^ + ' < 5 > 

Comparing Eqs. (4) and (5) we can identify L\\ = Lff, L\\ = L{\Z\/v z ni/vin z and 
L \ 2 = L \2 + L\\Zv z n t /v t n z . In Eqs. (4) and (5), c is the speed of light in vacuum and 
v a = ^2TJ m a is the thermal speed of species a (note that we employ Gaussian units). 
To facilitate comparisons between different models, we will perform a study at E r — 0 in 
this section but we will also use a more realistic finite value. The quantities B 0 , G, I 
and t in Eqs. (4)-(5) stem from the magnetic geometry specified in Boozer coordinates 0 
and ( in which 

b = k (V>, e , 0 VV> + / w Vd + G (^) VC. (6) 

B 0 is the (0,0) Fourier mode amplitude of B (9X), cI/2 is the toroidal current inside 
the studied flux surface, cGj 2 is the poloidal current outside the flux surface and t is the 
rotational transform [3]. In a typical stcllarator, |/| -C |G| and G ~ B 0 R , where R is the 
major radius of the device. Furthermore, in Boozer coordinates 



We study a W7-X vacuum configuration with i = — 1 at the plasma edge. For this 
geometry, the normalization factors of Eq. (5) are defined such that they are all positive 
except i and I (which have the same sign), and a positive flux is directed outwards whereas 
the density and temperature gradients in Eq. (3) are negative in the usual situation. The 
impurity transport coefficients, Lff, L\\ and Lj 2 , are obtained with SFINCS, solving 
coupled linear drift-kinetic equations for each species with three different right-hand- 
sides. Similarly to Ref. [14] we calculate the transport coefficients in terms of a normalized 
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collisionality for the impurities 


where 
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The collision operator implemented in SFINCS is discussed in Appendix A. In W7-X, the 
normalized collisionality can be expected to range from u' z ~ 0.01 in the core of a high- 
temperature, low-density, low-Z e g plasma to u' z ~ 10 at the edge of a low-temperature, 
high-density, high- Z e s plasma. In our calculations we will include unrealistically high 
values of i/ z to be able to compare our results with analytic theory which is only available 
at high collisionality (see Appendix C). It is interesting to note that collisions among the 
impurity ions themselves are more important than collisions with the bulk ions even if Z e s 
is not far above unity. This is because the ratio of the pitch-angle-scattering frequency 
between impurities and bulk ions i/f> and the pitch-angle-scattering frequency between 
impurity ions z/ff scales as ~ (Z 2 n z y / rn z j [16]. Thus, as soon as Z eS — 1 

significantly exceeds \J m*/ m z , z-z collisions are more important than z-i collisions. For 
a hydrogen plasma with C 6+ impurities, \/m % /rn z = 0.29. 

We also define a normalized electric field similar to Ref. [14], 

cG dd> 
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lv z Bq dxjj 

This is the electric field normalized by the so-called resonant electric field. In axisymmetry 
it corresponds to the poloidal Mach number. 


Simulation results at E r = 0 

In Figs. 1 and 2 the carbon transport coefficients as functions of v' z for the W7-X 
standard configuration geometry at r/a — 0.88 are shown, r is the effective radius related 
to the flux label through i/>n = 4’/4’a — ( r / a ) 2 > where a = 0.51 m is the outermost effective 
minor radius and V ; o, — V’ (Vfv = 1). At this radius the magnetic geometry parameters are 
B 0 = 3.1 T, G = 17.9 Tm, / = —6.5 x 10 -7 Tm and l = —0.93. The minimum resolution 
used in the SFINCS runs is N e = 17, = 49 grid points in the poloidal and toroidal 
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direction (per identical segment of the stellarator, where W7-X has a five-fold symmetry in 
the toroidal direction), N x = 5 grid points in energy (x = v/v a with v a being the thermal 
speed of the species a) and = 24 Legendre polynomials to represent the distribution 
function (here £ = vn/v), and Ni = 4 Legendre polynomials to represent the Rosenbluth 
potentials. Note however that the required resolution depends on the collisionality regime. 
At low collisionality and both typically need to be larger than 100, because of the 
presence of an internal boundary layer between the trapped-passing boundary. At high 
collisionality instead N x typically has to be increased. 

Results are presented for SFINCS simulations with full linearized Fokker-Planck- 
Landau collisions, and for Z e g = 1.05 at A* = 0 (see Fig. 1) and Z eS = 2.0 at 77* = 0 (see 
Fig. 2). Note that plots (b)-(d) sometimes use a double-logarithmic vertical scale, since 
these transport coefficients can have either sign. The results are compared to simulations 
with the DKES (Drift Kinetic Equation Solver) code [17, 18]. In contrast to SFINCS, 
DKES has 3 rather than 4 coupled phase-space coordinates, because energy coupling is 
neglected when solving the drift-kinetic equation. DKES employs pitch-angle scattering 
collisions, but momentum correction can be applied afterwards [13]. Here we use the 
momentum correction approach described in Ref. [11], This approach is appropriate for 
relatively collisionless cases, because the energy scattering part of the collision operator 
only contains the first order Legendre components of the distribution functions. Moreover, 
in the parallel particle and heat flows the approach neglects terms corresponding to the 
Pfirsch-Schliiter form of the poloidal E x B-drift. Consequently we will only show results 
at low collisionality (v' z < 1) using this momentum correction technique. We note that 
DKES employs different effective particle trajectories than SFINCS (compare Eqs. (17) 
and (18) in [14]), but the difference vanishes when E* = 0. In the short-mean-free-path 
limit, v' z > 1, the impurity transport coefficients can be computed analytically in terms 
of the parallel current. The details are given in Appendix C and based on the theory 
presented in Ref. [19]. We note that SFINCS retains E x B-precession when solving 
the drift-kinetic equation for f a \ [14], although according to the formal ordering it should 
appear first at next order. Since E x B-precession is not included in the analytical 
high-collisionality calculation in Ref. [19], we only compare SFINCS results to this liigli- 
collisionality theory at E* = 0. The high-collisionality asymptotes for the W7-X case are 
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shown in Figs. 12 and 13 in Appendix C. These theoretical limits conform well with the 
SFINCS computations in the appropriate limit. 


♦ Z e{{ = 1.05, E*=0, SFINCS Fokker-Planck-Landau x—x Z eff =1.05, E„=0, DKES + momentum correction 





FIG. 1: Carbon (Z = 6) transport coefficients Lff (a), Ll\ (b), Lf 2 (c) and the temperature 
gradient coefficient Zf 2 — 3 + Lf\j /2 (d) as functions of normalized collision frequency 

n' for a W7-X geometry at A* = 0 and with Z e g = 1.05. SFINCS computations with the 
Fokker-Planck-Landau collision operator ( —♦— ) are compared to DKES computations with 
momentum correction applied afterwards (— X— ). Note the double-logarithmic scale in (b)- 
(c). 


For both L\\ and Lf 2 the momentum correction technique captures well the intrinsic 
momentum conservation of the full linearized Fokker-Planck-Landau operator. Interest¬ 
ingly, the DKES + momentum correction curves fail to predict the sign change and 
1/n'-scaling of L\\ observed at low collisionality in the SFINCS curves with Fokker- 
Planck-Landau collisions of both Figs. 1 (b) and 2 (b), i.e. at E * = 0. 

With Fokker-Planck-Landau collisions at E * = 0, all impurity transport coefficients 
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♦ Z eff =2.0, E*= 0, SFINCS Fokker-Planck-Landau x—* Z eff =2.0, E*= 0, DKES + momentum correction 



FIG. 2: Carbon (Z = 6) transport coefficients Lff (a), Lf\ (b), L| 2 (c) and the temperature 
gradient coefficient Z| 2 — 3 + Lf\ \ /2 (d) as functions of normalized collision frequency 

v' z for a W7-X geometry at E* = 0 and with Z e s = 2.0. SFINCS computations with the 
Fokker-Planck-Landau collision operator ( —♦— ) are compared to DKES computations with 
momentum correction applied afterwards ( — X — ). Note the double-logarithmic scale in (b)- 
(d). 

show trends of l/^'-transport at low collisionality and are proportional to z/ at high col- 
lisionality. Earlier work had assumed that the inter-species transport coefficients should 
be negligible compared to the self-species coefficients at low collisionality since the species 
interact via collisions [2], This is not in agreement with what we find in our SFINCS cal¬ 
culations at E* = 0 as shown in Figs. 1 and 2, where the Fokker-Planck-Landau curves 
of L\\ also exhibit l/z/'-transport at low collisionality. This phenomenon is explained as 
follows. Since L^\ is found by setting all thermodynamic forces to zero except the main 
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ion gradient, the relevant impurity kinetic equation is 

V\\ V||/ Z 1 = C zz [fzl, fMz] + C zz [fMz, fzl] + C z i [fzl, fMi] + C z i [/ Mz, fil] , (H) 

where C a b [fa, ft] is the Fokker-Planck-Landau collision operator for species a with dis¬ 
tribution f a colliding with species b of distribution /),. At low collisionality the particle 
transport is carried by the trapped particles. If we perform bounce-averaging we annihi¬ 
late the streaming term in Eq. (11) and obtain 

0 = C zz [fzl, f Mz\ + Czz [fMz, fzl] + C zi [f z i, f M i] + C zi [f Mz, fil], (12) 

where the bounce-average is denoted by overhead bar. In Eq. (12) every term contains a 
factor v' z which can be divided away. The equation is a linear inhomogeneous equation 
f° r fzi, where the inhomogeneous drive term is C zi [fMz, fa] containing fn and which is 
non-zero because of the cZrz*/drive in the main ion kinetic equation. Because of the 
1/v- regime of a stellarator in the absence of E *, we expect that fn scales like 1 jv' z and 
it could thus be expected from Eq. (12) that also f zi scales as 1 /v' z giving rise to the 
behavior in L\\ we observe at low collisionality. Note that the l/i/-part of fn is even in 
un, and since U|| parity is preserved by the field term of the linearized collision operator, 
then the even part of the C Z i field term is required to couple this drive to the impurities. 
Hence, the l/z/-scaling of L\\ will be missed in any numerical or analytic calculation in 
which the terms that are even in U|| are neglected in the field term of the collision operator. 
Momentum conservation, which is associated with the odd part of the field term, is not 
sufficient. For this reason, the momentum-corrected DKES results in Figs. 1 (b) and 
2 (b) obtain the wrong scaling with v' z and wrong sign at low collisionality. It is possible 
to test this hypothesis using a modified form of the impurity-ion Fokker-Planck-Landau 
collision operator in SFINCS, by selectively turning off the field term for even Legendre 
modes. The results of these calculations for Z e s — 1 .05 at E* = 0 are shown in Fig. 3, and 
clearly the l/i/-behavior of Lf\ disappears when the even Legendre polynomials in the 
field term of C Z i have been suppressed. Physically, if the bulk ions have a radial density 
gradient, then fn carries an anisotropic pressure and will be rich in particles drifting 
outwards and poor in particles drifting inwards. The term C zi [$m z , fa] will try to create 
a similar anisotropy in f z and thus cause radial impurity transport. 
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4-^ Z e!! = 1.05, E,=0, SFINCS Fokker-Planck-Landau 



FIG. 3: Comparison of the Ll\ coefficient as function of normalized collision frequency u' z 
calculated with SFINCS using full linearized Fokker-Planck-Landau collisions and when the 
even Legendre polynomials in the field term of the impurity-ion collision operator C z i have been 
suppressed, for a W7-X geometry at If* = 0: Full linearized Fokker-Planck-Landau operator and 
Z e g = 1.05 ( —♦— ), even Legendre polynomials in field term of C z i suppressed and Z c g = 1.05 
(-*-)■ 

A negative L\\ at all collisionalities is not surprising, since it merely tells us that 
a negative impurity density gradient will drive the impurities outwards. This is nec¬ 
essary and follows from the entropy law. Moreover, as shown in Appendix C, in the 
high-collisionality limit Lf x = — ZL\\ . (We note that this is also approximately valid for 
the high-collisionality SFINCS Fokker-Planck-Landau results at A* = —0.74 in Fig. 4, 
although SFINCS retains E x B-precession which is formally excluded in the usual drift 
ordering.) This implies that the main ion density gradient will drive an impurity accu¬ 
mulation, which is stronger the higher the impurity charge. Furthermore, in accordance 
with the discussion for tokamaks in Ref. [16], in the Phrsch-Schliiter regime the impurity 
transport is primarily driven by the bulk ion gradients. In the low-collisionality banana 
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regime, however, L\\ can be smaller in size than L\\ and also be negative. This indicates 
that the bulk ion density gradient could mitigate an impurity accumulation in a hot re¬ 
actor plasma. Note that since L\\ = —ZL\\ in the Pfirsch-Schliiter limit, and because of 
our definition of the thermodynamic forces and the transport coefficients in Eqs. (3) and 
(5), the ( Ze/T z ) (dQ/di/j) term in A zl is canceled by the (e/7)) term in A tl in 

Eq. (5). The reason why the radial electric field has no effect on the impurity transport in 
this limit is that the transport is dominated by impurity-ion friction, and transport from 
friction is intrinsically ambipolar [19]. However, it is important to emphasize that this is 
a result in the usual drift ordering where the E x B-drift velocity is ordered Ve ~ p*Vi and 
E x B-precession is formally excluded. Recall that SFINCS (as well as other stellarator 
neoclassical calculations [17, 18]) retains E x B-precession in the drift-kinetic equation 
for f a i, and it is thus possible that a finite E * yields a different Pfirsch-Schliiter limit 
than E* = 0 in SFINCS calculations. As noted in Ref. [10], poloidal E x B-precession 
is traditionally ignored in the Pfirsch-Schliiter regime but can become relevant when the 
product of the radial electric field and collisionality is sufficiently large. 

The L\ 2 ~ 3 + Ll\J /2 coefficient illustrated in Figs. 1-2 (d) is the temperature 

gradient coefficient (recall that we assume T) = T z ), i.e. the coefficient in front of dlnT/di/j 
when substituting Eq. (3) into Eq. (5). This coefficient corresponds to — ( D 2 Z + D 2 ) in 
Eq. (1), adjusted with a normalization factor. A negative value indicates temperature 
screening which is found in the low-collisionality regime. Figure 2 (d) shows that in 
the high-collisionality regime SFINCS finds that the temperature gradient drives the 
impurities inwards, with a strength increasing with collision frequency when Z e g = 2.0 
and E* = 0, which is in agreement with the analytical Pfirsch-Schliiter asymptotes (see 
Appendix C). For Z e s = 1.05 and E* = 0 the SFINCS Fokker-Planck-Landau calcula¬ 
tions show a temperature screening also in the high-collisionality regime, as illustrated 
in Fig. 1 (d). Our SFINCS results are again in agreement with the analytical Pfirsch- 
Schliiter asymptotes, and it is intriguing to see that whether we find a temperature screen¬ 
ing in the Pfirsch-Schliiter limit or not can, in fact, depend on the impurity content. It 
should be emphasized that these extremely high collisionalities are irrelevant in practice 
for at least two reasons. Firstly, the observed transport is usually turbulent in plasmas 
that are cold enough to be in the Pfirsch-Schliiter regime; and secondly, the neoclassical 
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transport is sensitive to the radial electric field. 

Since Fff, L*\ and L z l2 — 3 + Lf, j /2 (i.e. the impurity density gradient, the 

ion density gradient and the temperature gradient coefficients respectively) are all nega¬ 
tive in the low-collisionality regime at F* = 0, one might think that both the main ion 
and impurity pressure gradients would be beneficial for avoiding neoclassical impurity 
accumulation in a hot reactor-like plasma if the profiles are peaked. However it should 
be remembered that the radial impurity transport will typically heavily depend on the 
ambipolar radial electric field which builds up to balance the particle fluxes and to yield 
a vanishing radial net current. 


Simulation results at finite E r 

To examine the role of the radial electric field, we perform SFINCS calculations 
with full linearized Fokker-Planck-Landau collisions and Z e g = 2.0 at F* = —0.74. This 
value of F* corresponds to E r = —20 kV/m if the temperature is T = 1 keV. The results 
are illustrated in Fig. 4 and when comparing the SFINCS curves to the DKES curves 
it should be recalled that the simulations employ different effective particle trajectories, 
which matters when F* ^ 0. To examine if the difference between the SFINCS results 
and the DKES results is mainly a consequence of the different collision operators or 
the different effective particle trajectories which are employed in the two tools, we have 
also included results from SFINCS calculations using “DKES particle trajectories”. As 
shown in Fig. 4 the effect of the different models for the particle trajectories is small 
for this particular case. However, note that as F* approaches unity a difference can in 
general be expected [14]. In contrast to the results for F* = 0, at F* = —0.74 there is 
no sign change in the SFINCS Fokker-Planck-Landau curves for L\\ at low collisionality, 
as seen in Fig. 4 (b), but there is still a difference to the DKES + momentum correction 
results. Moreover, the 1/z/'-transport at low collisionality and i/-proportionality at high 
collisionality seen in Figs. 1 and 2 disappear for the finite F*. The temperature coefficient 
shown in Fig. 4 (d) indicates a temperature screening at all plotted collisionalities, but 
the effect is weak at low collisionality. In the low-collisionality regime we also find a 
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negative impurity density gradient coefficient, but a small positive ion density gradient 
coefficient. The magnitude of L\\ is significantly smaller than the magnitude of Lff, which 
implies that if we substitute a negative E r into Eqs. (3) and (5) the electric field will cause 
impurity accumulation as expected. 

♦—♦ Z e(l = 2.0, E*=- 0.74, SFINCS Fokker-Planck-Landau 
*—* Z eE = 2.0, E t =-0.74, SFINCS Fokker-Planck-Landau, DKES trajectories 



FIG. 4: Carbon (Z = 6) transport coefficients Lff (a), Lf\ (b), L\ 2 (c) and the temperature 
gradient coefficient Zf 2 — 3 + Lf\j /2 (d) as functions of normalized collision frequency u' z 

for a W7-X geometry at E* = —0.74 and with Z e g = 2.0. SFINCS computations with the 
Fokker-Planck-Landau collision operator ( —♦— ) are compared to DKES computations with 
momentum correction applied afterwards ( — X — ). Also shown are SFINCS computations 
using a different model for the particle trajectories, referred to as “DKES particle trajectories” 
and described in Ref. [14] ( —★— ). Note the double-logarithmic scale in (c). 

As a final remark of this section, we note that we can also use SFINCS to calculate 
the main ion transport coefficients in Eq. (4). We find that L\\ = Lf\ (both for finite 
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and vanishing A*), and thus the Onsager symmetry described in the beginning of the 
section is fulfilled. 


III. IMPURITY DENSITY PEAKING AND BOOTSTRAP CURRENT 

In a 3D device, a bootstrap current arises due to similar reasons as in a tokamak, 
but the size of it is typically substantially smaller than the Ohmic current in a tokamak 
[1, 2, 20]. The bootstrap current is a consequence of the trapped particle orbits, and is 
generally larger at low collisionality than at high collisionality [16, 21], Since the bootstrap 
current adds pressure-dependence to the magnetic equilibrium, it can be desirable in 
stellarators to minimize the bootstrap current so the magnetic field remains optimized 
over a range of plasma pressure. W7-X has been optimized for a small bootstrap current. 
A net toroidal current changes the value of l at the boundary, which can be detrimental 
for proper island divertor operation [22-24], It is consequently important to be able to 
make realistic predictions of the bootstrap current when designing a stellarator. 

In this section we use SFINCS to investigate how the presence of impurities affects 
the bootstrap current in a non-axisymmetric plasma, and also how the neoclassical impu¬ 
rity dynamics are affected by the impurity content through the plasma effective charge. 
Again we study a hydrogen plasma with a single carbon impurity species present (this 
time the electrons are included in our simulations), and use the W7-X standard magnetic 
configuration. In contrast to Sec. II, the configuration here corresponds to a central elec¬ 
tron cyclotron resonance heating profile and an average /3 of 2.9 %, /3 = 8imT / B 2 being 
the ratio of plasma pressure to magnetic pressure. The temperature profiles have been 
calculated according to the procedure in Ref. [25], in which the electron density profile 
is assumed on the basis of other experiments. The density and temperature profiles are 
shown in Fig. 5. 

To study the radial impurity transport, we calculate the neoclassical zero-flux impu¬ 
rity density gradient (also referred to as the impurity peaking factor) defined as ci/L nz for 
which the flux-surface-averaged impurity flux vanishes, (T 2 ■ V?/’) = 0. This happens when 
the convective part of the transport is balanced by the diffusive part, corresponding to V z 
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r/a 


r/a 


FIG. 5: The plasma profiles for a W7-X standard magnetic configuration corresponding to 
central electron cyclotron resonance heating profile and an average /3 of 2.9 %: (a) assumed 
electron density (-); (b) predicted electron ( ) and ion temperatures ( - ). 

and D\ z in Eq. (2). L nz = — [d\nn z /dr]~ l is the impurity density gradient scale length. 
(Note that Lt = — [d In T/dr]^ 1 is kept fixed when calculating the zero-flux impurity den¬ 
sity gradient.) To calculate the impurity peaking factor we must also find the ambipolar 
radial electric field simultaneously, i.e. the E r for which the radial net current vanishes, 
Z^e (T& • Vr/’) = 0. As earlier mentioned, in a typical ion root scenario (negative E r ) 

b 

where the impurity concentration is small compared to the main species, the ambipolar 
electric field arises to bring the main ion particle transport down to the electron level, 
and so the ambipolar electric field is approximately the field for which (T* ■ V^>) ~ 0. For 
a single impurity species plasma, with the impurity species in trace contents, it is always 
possible to find values of driz/d'ip and E r such that the radial impurity flux and the radial 
current vanish simultaneously. We note however, that if Z e g = Z (i.e. a plasma consist¬ 
ing of only electrons and the impurity species) ambipolarity requires the impurity flux to 
balance the outward electron flux, and it is not possible to find a neoclassical impurity 
peaking factor. Thus, there must be a critical value of Z e g between 1 and Z above which 
the radial impurity flux and the radial current cannot vanish simultaneously. (This is the 
reason why we later in this section are not able to present the impurity peaking factor at 
large values of Z eS for one of the cases, i.e. no solution exists.) 

We analyze two radial locations, r/a = 0.2 and r/a = 0.8. At r/a = 0.2 the 
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parameters are n e = 0.991 x 10 20 m -3 , T e = 5.44 keV, T % = 3.60 keV, dn e /dr = —0.209 x 
10 20 m -4 , dT e /dr = —13.5 keV/m and d:T { /dr = —0.269 keV/m. At r/a = 0.8 they 
are n e = 0.758 x 10 2 ° m -3 , T e = 1.74 keV, Tj = 1.82 keV, dnjdr = —6.93 x 10 20 m -4 , 
dT e /dr = —13.1 keV/m and dTi/dr = —13.1keV/m. Z e g is varied by varying the impurity 
density (and main ion density accordingly to fulfill quasi-neutrality at hxed n e ), and this 
implies that the normalized collision frequency, defined in Eq. (8), satisfies u z < 0.25 
for r/a = 0.2 and v' z < 0.70 for r/a = 0.8 respectively. It is difficult to anticipate on 
what time scales the impurity species will reach a steady-state (and its corresponding 
zero-flux impurity density gradient), and also how the main species density gradients will 
adapt to keep radial quasi-neutrality. We therefore carry out two different scans in Z e ff. 
In the first scan we find the impurity peaking factor, and the main ion density gradient 
is modified along with the impurity density gradient to preserve radial quasi-neutrality, 
while the electron density gradient is kept fixed. (The corresponding curves are solid and 
labeled “zero-flux carbon gradient” in Fig. 6.) In the second scan we instead keep the 
radial density gradient scale lengths fixed and equal, L ne = L ni = L nz , which is equivalent 
to the condition dZ eS /dr = 0. (The corresponding curves are dashed and labeled u Z eS 
independent of r” in Fig. 6.) 

Figure 6 shows our SFINCS results for the carbon density gradient, ambipolar radial 
electric held, bootstrap current density, and carbon how as functions of Z eS for the W7-X 
geometry. The scans are performed with full linearized Fokker-Planck-Landau collisions. 
Comparing the approach when we fold the carbon peaking factor to the approach when we 
keep the density gradient scale lengths hxed, we see that qualitatively a similar ambipolar 
electric held and bootstrap current density are obtained (compare the solid and dashed 
curves in Figs. 6 (b) and (c)). Unsurprisingly the deviation between the two approaches 
increases as Z e g becomes larger, mainly because the difference in main ion density gradient 
also becomes larger. A similar argument, but regarding the impurity density gradient, 
is likely the reason why there is a difference in parallel carbon how between the two 
approaches (see Fig. 6 (d)). Note that particularly the core carbon how (both the direction 
and the magnitude) at r/a = 0.2 is sensitive to Z e g. Thus the how could make a sensitive 
diagnostic test of neoclassical physics. 

Irrespective of which model for the density gradients we use, some main features 


20 


are clearly observed. Firstly, at r/a = 0.2 the plasma is in the electron root regime 
and at r/a = 0.8 in the ion root regime, since the ambipolar radial electric field has the 
opposite sign as shown in Fig. 6 (b). Furthermore, in accordance with what was found 
in Ref. [26] the ambipolar electric field is reduced as the impurity content is increased. 
This also implies that the impurity profile, which is hollow at r/a = 0.2 and peaked 
at r/a = 0.8, is flattened as illustrated in Fig. 6 (a). The bootstrap current density is 
also significantly reduced with increased impurity content, although no sign change is 
observed (see Fig. 6 (c)). This reduction is not surprising, since the electron-ion friction 
increases with Z e s- The results show that changes in Z e g on the order of AZ e g ~ 1 can 
lead to changes in the bootstrap current density on the order of Ajbs > 20kA/m 2 . If such 
changes occurred across the entire plasma cross-section, the total current could change by 
A I hs > 10 kA. This is illustrated in Fig. 7, where the integrated total bootstrap current 
is shown as a function of Z e g when Z e g is kept constant throughout the radial domain. A 
change in the plasma current of 10 kA could modify the value of the boundary -1 enough to 
cause measurable changes in the island divertor strike point locations [22-24], indicating 
that the full ion composition should be considered when performing bootstrap current 
calculations. 


IV. DISCUSSION AND CONCLUSIONS 

In this work we have used a continuum drift-kinetic solver, the SFINCS code, to 
calculate neoclassical impurity transport coefficients (defined by Eqs. (3) and (5)) in a non- 
axisymmetric magnetic equilibrium. Particularly, we studied carbon transport close to the 
plasma edge in a W7-X hydrogen plasma for two different levels of the impurity content 
corresponding to Z eS = 1.05 and Z eS = 2.0 respectively, at vanishing radial electric field. 
For Z e ff = 2.0 we also studied the carbon transport at E * = —0.74 (corresponding to 
E r = —20 kV/m if the temperature is T = 1 keV). 

We compare SFINCS computations with full linearized Fokker-Planck-Landau colli¬ 
sions to computations with the DKES code (which employs pitch-angle scattering) where 
momentum correction is applied afterwards. We find that the impurity density gradient 
and temperature gradient coefficients are well reproduced by the numerical model with 
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■■r/a = 0.2, Z e ff independent of r, Fokker-Planck-Landau 
•r/a = 0.8, independent of r, Fokker-Planck-Landau 


-r/a = 0.2, zero-flux carbon gradient, Fokker-Planck-Landau 
-r/a = 0.8, zero-flux carbon gradient, Fokker-Planck-Landau 


Impurity gradient 
a/L nz 


E r [kV/m] 


Bootstrap current 
- O' • B) /B 0 [kA/m 2 ] 




Carbon flow 
(Vc ■ B) /Bo [km/s] 



FIG. 6: Carbon (Z = 6) density gradient (a), ambipolar radial electric field (b), bootstrap 
current density (c), and carbon flow (d) as functions of the plasma effective charge for a W7-X 
geometry. SFINCS computations with Fokker-Planck-Landau collisions at two different radii 
are compared, using two different models for the ion density gradients, one when the carbon 
peaking factor is found (solid lines labeled “zero-flux carbon gradient”) and another when Z e Q is 
kept independent of the radial location (dashed lines labeled “Z e g independent of r”): r/a = 0.2 
(+), r/a = 0.8 (x). 

pitch-angle scattering and momentum correction, but to correctly determine the inter¬ 
species ion density gradient coefficient it is sometimes not sufficient to merely account for 
momentum conservation. 

The impurity transport coefficients show trends of l/z/-transport at low collision- 
ality and are proportional to v' z at high collisionality if E * = 0. From earlier work 
l/z/'-transport is not expected for the ion density gradient coefficient. We show that, 
by suppressing the even Legendre polynomials in the field term of the impurity-ion colli¬ 
sion operator in SFINCS the l/i/-transport for this coefficient disappears. Earlier work 
has often approximated the field-particle part of the collision operator by a momentum- 
conserving term, which has the wrong vy parity to couple the (even) l/z/'-part of fn to 
the impurities and this is likely the reason why 1 /z/'-transport for this cross-species trans- 
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Bootstrap current density Total bootstrap current 



FIG. 7: Radial profiles of the bootstrap current density from SFINCS calculations with the 
full linearized Fokker-Planck-Landau collision operator at different Z e s (keeping Z e g constant 
throughout the radial domain) and the integrated total bootstrap current as a function of Z e g. 

port coefficient has not been observed. If we introduce a finite E* in onr calculations, we 
instead find trends of -y/fZ-transport at low collisionality. Not surprisingly, we find that 
the impurity density gradient coefficient is negative at all collisionalities which merely 
implies that a standard negative gradient drives the impurities outwards. At high colli¬ 
sionality however, the impurity transport is dominated by the bulk ion density gradient 
(by a factor Z) which drives the impurities inwards. At low collisionality we find that all 
transport coefficients are negative when E * = 0, indicating an impurity screening. This 
could be beneficial from a reactor point of view, since the hot (almost collisionless) core 
could avoid impurity accumulation. However, when we introduce a negative ambipolar 
radial electric held (ion root regime) it drives impurity accumulation as expected. Inter¬ 
estingly we fold a temperature screening in the low-collisionality regime which persists up 
to a relatively high collisionality for Z c ft — 2 .0, and is maintained at all collisionalities for 
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Z e s — 1 .05. This is also the case for the calculations with a finite E *. In the high colli- 
sionality limit, an analytic prediction is available for the transport coefficients at E t = 0, 
and the SFINCS calculations conform well with these predictions. 

Moreover, we have used SFINCS to investigate how the impurity content affects the 
neoclassical impurity dynamics and the bootstrap current in a W7-X plasma. We find that 
an increased impurity content, implying a higher plasma effective charge, tends to flatten 
the impurity profile (determined from the condition of zero impurity particle flux) both 
close to the core (r/a = 0.2) where it is hollow and close to the edge (r/a = 0.8) where 
it is peaked. This trend is attributed to the reduction of the ambipolar radial electric 
held with increasing Z e s- The bootstrap current is also reduced with increasing impurity 
content, which is expected since the electron-ion friction increases with Z e s. Importantly, 
we find that the change in bootstrap current can be larger than 10 kA for a change in Z e s 
of O (1). A change of this size could be large enough to cause a deviation in the divertor 
strike point locations. This emphasizes the importance of performing bootstrap current 
calculations with a realistic ion composition. 
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A. COLLISION OPERATOR IN SFINCS 


The total collision operator for species a is a sum of linearized collision operators with 
each species: C a [f a ] = ^ C l ab [f a , f b ] , where C l ab [f a , f b ] = C ab [/„ i, f Mb ] + C ab [f Ma , hi] 

b 

and C ab [f a , f b ] is the Fokker-Planck-Landau collision operator. C ab [f a i, fMb] is referred 
to as the test particle part and C ab [/m q , f b i\ is the field particle part [9, 16, 27]. The 
linearized collision operator may be written as C l ab = C^ b + C^ b + C^ b , where C^ b + C^ b to¬ 
gether represent the test particle part and the field particle part. The Lorentz part is 


Cah = VDab (v) C [f a \ with V Dah ( V) = U ab 


<t> (x b ) - T (x b ) 


being the deflection frequency and 

1 d d 2 f x 

C = (l — £ 2 ) — the Lorentz operator. In this context (f) (x) = —j= / exp (— y 2 ) dy 

2 at, ot, y tt J o' 


is the error function, T (x) = 


(f) (x) — x<j>' (x) . 


2x 2 


is the Chandrasekhar function, v ab = 


2 1 / 2 7t n b e 2 a el In A/ and x a = v/v a with v a = Tj m a being the thermal speed 

of the species a. The energy scattering contribution is 


eft = mn> 


d 2 f, 


al 


2 dv 2 


Xu 


1 _ y-f^l 
m b J dv 


+ ^DabV-^- + 4:'K—0 ab — f Mb fal, (Al) 

dv n b m b 


V(x b ) 


where v\\ ab = 2v ab --— . We can write the field term as 
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Ch = Ch + C^ + C, 
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(A2) 


where 


/VI _ ^ab r 

L 'ab ~ JMa 

n b 



dH b i 
dv 


2v a H b \ 


r G - 



n b 


2v 2 d 2 G b i 
v a dv 2 


(A3) 

(A4) 


and 

Cab = — vlfMa^n—hi- (A5) 

n b m b 

The functions G b i and H b 1 are the perturbed Rosenbluth potentials [9] defined by V 2 H b i = 
-47 xfbi and V 2 G W = 2 H bl . 

The speed discretization in SFINCS is based on a spectral collocation scheme de¬ 
scribed in [27]: a function f (x) is stored at grid points x 3 which are the zeros of a 
polynomial, where the polynomial is taken from the set (x) (with n > 0) obeying the 
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orthogonality relation 

dxx k exp (-x 2 ) M k (x) ( x ) = 8 n , m A k n . (A6) 

Here S n ^ m is the Kronecker delta, A k represents some normalization and k is any number 
greater than —1, but from experience the choice k = 0 is typically good. Note that x here 
denotes the speed normalized to thermal speed, e.g. the distribution function for species 
a is f a (x a ). We can alternatively represent / in a modal discretization by the vector of 
numbers F k in 

/ (x) = KK M ex P (—x 2 ), (at) 

n 

where 

i r°° 

dxx“Ml(x)f(x). (AS) 

JO 

In terms of the Gaussian integration weights Wj associated with the grid Xj satisfying 

OO 

dx y(x) ~ Wjy[xj ), there is thus a linear transformation Y from the collocation 

3 

to the modal discretization, with Y m ^ n = w n x k M^(x n ) / A^. 

At each speed grid point, the pitch-angle dependence of the distribution function is 
decomposed in Legendre polynomial modes: 



fal ( x a, £) = ^ p i (0 fal,l M • (A9) 

l 

In the Legendre modal representation, C r ah becomes diagonal. As in Ref. [27], Cff b can be 
represented using the pseudospectral differentiation matrix associated with the polyno¬ 
mials M k (x), and Cjf b can be represented using the interpolation matrix associated with 
the M k (x). In evaluating Cff and Cf h , however, we depart from the method in [27]. First 
an expansion analogous to Eq. (A9) is made for the perturbed potentials in terms of their 
Legendre modes H b \j (x b ) and G b i t i (x b ), and we find 


_d_ 2 dH bl>l 
dxb Xb dxb 


1(1 + 1) H bl) i = -47 rv 2 f b ij, 


d^dGbu _ /(/ _ |} Gbii 2 v 2 H bhl . 

OX 5 OX 5 

Solving Eqs. (AlO)-(All) with a Green’s function approach, we obtain 


Hbl,l ( x b) — 


47T 


(2(+ l) U! +1 Jo 


r x b 
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and 


Gbl,l (x b ) — 


47T 


(4/2-1) 
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dz z l+i fbi,i (z) - 
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dzz l+4 f bli i(z) + -j-f I dzz l+2 f bli i{z ) 




(A13) 


as integrals of / 61>i , as in Eqs. (40) and (45) in Ref. [9]. To find and we need 
dH b i i/dx b and d 2 G b i y i/dxl which are computed by analytically differentiating Eqs. (A12)- 
(A13). We evaluate Eqs. (A12)-(A13) and their derivatives replacing f blti by each of the 
polynomials M%, using integration endpoints x b corresponding to each speed grid point 
for species a normalized to v b . For each Legendre mode, the results for N x polynomials 
and N x evaluation points yield a N x x N x matrix which we denote by R. Thus, the 
map from distribution for species b (on the speed collocation grid points for species b ) 
to perturbed Rosenbluth potentials (on the speed collocation grid points for species a) 
is given by the matrix product RY. The computational expense of these integrations is 
negligible compared to solving the main linear system of discretized kinetic equations. For 
most circumstances, the method of evaluating C^ b and C^ b described here gives identical 
results (to 2 or more decimal places) to the method in [27]; however we find the method 
here to yield better convergence at extremely high collisionality. 


B. IMPURITY TRANSPORT COEFFICIENTS FOR PITCH-ANGLE 

SCATTERING MODELS 

In Figs. 8, 9 and 10 we compare the SFINCS Fokker-Planck-Landau computations 
of the impurity transport coefficients in Sec. II to SFINCS computations with pitch- 
angle scattering and DKES computations (no momentum correction applied afterwards). 
Moreover, in Fig. 11 we compare SFINCS pitch-angle scattering computations using 
“DKES particle trajectories” to DKES computations. It is reassuring to see that when 
the two different numeric tools use the same collision operator and the same effective 
particle trajectories, they yield practically the same results. We also see that for this par¬ 
ticular case the resulting difference from using different models for the particle trajectories 
is small. 
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Figures 8-10 (a) and (c) show that at low collisionality, momentum conservation is 
unimportant for L\\ and L\ 2 . This finding is consistent with the results of [14] where 
a single ion species was analyzed, and is explained as follows. In the low-collisionality 
1/z/-regime the radial transport is connected to pitch-angle scattering of helically trapped 
particles, and the dominant physics is captured by the pitch-angle scattering approxima¬ 
tion. If a radial electric held is present this is also true for the y/v^- regime. The effect of 
the collisions is mainly to scatter particles across the trapped-passing boundary in velocity 
space. 

In the high-collisionality regime, the difference in L\\ is small between the 
momentum-conserving linearized Fokker-Planck-Landau operator and the pitch-angle 
scattering operator at low Z c g, whereas at Z e g = 2.0 and E * = 0, the pitch-angle¬ 
scattering result is a factor of ~ 5 larger than the Fokker-Planck-Landau result. However, 
at Z e g = 2.0 and E * = —0.74 the difference is smaller. In contrast, for L\ 2 the sign of the 
coefficient depends crucially on which collision operator is used for both values of Z e s- 
This is also verified by the results in Appendix C. For both L\\ and Lf 2 the DKES curves 
conform reasonably well with the SFINCS pitch-angle scattering curves. 

Furthermore, we see that for pitch-angle scattering the ion density gradient coef¬ 
ficient disappears in all collisionality regimes: L\\ = 0 (and also L z ]2 = 0), which is 
confirmed by the results presented in Appendix C and can be understood as follows. In 
the absence of a momentum-conserving term, the impurities only feel collisions with a sta¬ 
tionary background. In the impurity drift-kinetic equation, there is no information about 
the density gradient of the main (bulk) ions, and consequently the radial impurity flux is 
independent thereof. If momentum is conserved in the collisions, however, the impurities 
are affected (through the collision operator) by the bulk ion flux along the magnetic held, 
which depends on the ion density gradient. 

Finally, Figs. 8-9 (d) show that at intermediate and high collisionality, the absence 
of momentum correction in the collision operator can lead to transport predictions in 
the wrong direction for the temperature gradient coefficient. In the SFINCS Fokker- 
Planck-Landau calculations a temperature screening is typically found (except at very 
high collisionality for the Z e g = 2.0 case), but pitch-angle scattering calculations predict 
an inward impurity drive. However for the calculations at finite E* (Fig. 10 (d)) both 
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collision models give practically the same result for the temperature gradient coefficient, 
and a screening is found at all collisionalities. 


♦ Z eff = 1.05, E t = 0, SFINCS Fokker-Planck-Landau 
•—• Z e{{ = 1.05, E*= 0, SFINCS Pitch-angle-scattering 




FIG. 8: Carbon (Z = 6) transport coefficients L\\ (a), (b), L\ 2 (c) and the temperature 

gradient coefficient L\ 2 — 3 /2 (d) as functions of normalized collision frequency n' 

for a W7-X geometry at hi* = 0 and with Z e g = 1.05. SFINCS computations for two different 
collision operators are compared: Fokker-Planck-Landau ( —♦— ) and pitch-angle scattering 
( —•— ). Also shown are results from DKES (pitch-angle scattering) ( —■— ). Note the double- 
logarithmic scale in (b)-(d). 


C. IMPURITY TRANSPORT COEFFICIENTS AT HIGH COLLISIONALITY 

In Ref. [19] analytic calculations for the impurity transport in the Pfirsch-Schliiter 
regime are presented, and from these we can derive expressions for L\\ and L\ 2 . The 
drift-kinetic equation for the first order (in p*) distribution function f a \ is solved by a 
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♦ Z e{! = 2.0, E*= 0, SFINCS Fokker-Planck-Landau 
•—• Z eS = 2.0, E„ = 0, SFINCS Pitch-angle-scattering 
■—■ Z ef[ =2.0, E x = 0, DKES 





FIG. 9: Carbon (Z = 6) transport coefficients Lff (a), Lf\ (b), Lf 2 (c) and the temperature 
gradient coefficient L\ 2 ~ 3 (jXff + L*\j /2 (d) as functions of normalized collision frequency v' z 
for a W7-X geometry at E* = 0 and with Z e g = 2.0. SFINCS computations for two different 
collision operators are compared: Fokker-Planck-Landau ( —♦— ) and pitch-angle scattering 
(— •— ). Also shown are results from DKES (pitch-angle scattering) ( —■— ). Note the double- 
logarithmic scale in (b)-(d). 

subsidiary expansion in the shortness of the mean free path, Aj = Xu/L <C 1, where 
Xu = Vi/u is the ion mean-free-path and L ~ V -1 is the plasma dimension. The lowest 
order solution ~ A* -1 ) is a shifted Maxwellian. The impurity flux is determined 

from a pressure anisotropy term and an impurity-ion friction term, whose relative sizes 

pressure anisotropy term Xf, . 

are - - --- ~ ——-. In the short-mean-tree-path limit, the pressure 

friction term L 2 Z 4 

anisotropy term can consequently be neglected in transport calculations. The friction 
term is intrinsically ambipolar, and the fluxes are independent of the radial electric field 
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♦ Z eff =2.0, E* = - 0.74, SFINCS Fokker-Planck-Landau 
•—• Z e[[ =2.0, E x =—0.74, SFINCS Pitch-angle-scattering 
■—■ Z eff =2.0, E,=- 0.74, DKES 





FIG. 10: Carbon (Z = 6) transport coefficients Zff (a), iff (b), Zf 2 (c) and the temperature 
gradient coefficient Zf 2 — 3 ^Zff + Zff ^ /2 (d) as functions of normalized collision frequency 
v' z for a W7-X geometry at E„ = —0.74 and with Z e g = 2.0. SFINCS computations for 
two different collision operators are compared: Fokker-Planck-Landau ( —♦— ) and pitch-angle 
scattering (— •— ). Also shown are results from DKES (pitch-angle scattering) ( —■— ). Note 
the double-logarithmic scale in (c) and that the pitch-angle scattering results in (b) vanish. 


in this limit (when the usual drift ordering Ve ~ is used and E x B-precession is 
formally excluded). 

The coefficients are straightforwardly obtained from Eq. (2) in Ref. [19], after ne¬ 
glecting the pressure anisotropy term. Note that Ref. [19] employs Si-units, and we need 
to transform the corresponding expressions into Gaussian units to match Eq. (5). The 
impurity coefficients depend on the geometry-dependent quantity u satisfying 


Villi = 


4(bx VZ) ■ VlnR, 


(A 14) 
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Z eE = 2.0, £, = -0.74, SFINCS Pitch-angle-scattering 

Z eff =2.0, £„ =—0.74, SFINCS Pitch-angle-scattering, DKES trajectories 



FIG. 11: Carbon (Z = 6) transport coefficients Lff (a), Lf\ (b), L\ 2 (c) and the temperature 
gradient coefficient Zj 2 — 3 + Lf \^ /2 (d) as functions of normalized collision frequency v' z 

for a W7-X geometry at £* = —0.74 and with Z e g = 2.0. SFINCS computations for pitch- 
angle scattering with two different models for the particle trajectories (described in Ref. [14]) 
are compared: “full particle trajectories” ( —•— ) and “DKES particle trajectories” ( ~ ). 

Also shown are results from DKES (pitch-angle scattering) ( —■— ). Note that Ljj in (b) vanish. 


where Vm = b- V is the gradient along the magnetic field, u is proportional to the parallel 
current divided by B. 

Expressions for the impurity transport coefficients in the high collisionality regime, 
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and under the assumption T z — T t — T, are summarized in Eqs. (A15) and (A16). 



(A 15) 


W = - <« 2 B 2 > < °. (A16) 


(A16) 


ao, (3 q and are coefficients of expanded in Sonine polynomials, with being 


the zeroth order term in an expansion in the shortness of the ion mean-free-path of the 
first order (in p*) ion distribution f t \. The details are given in Ref. [19], and here we have 
calculated the coefficients for Z = 6 and m z /rrii = 11.924. For Z e s = 1.05 we obtain 
a 0 = —24.689, a n = 1.076, /3 0 = —2.689, = 2.857, and for Z eS = 2.0 we obtain 

a 0 = -1.741, a h = 0.574, (3 0 = -1.435, Ai = 1-812. 

From a calculation similar to Ref. [19] but only including pitch-angle scattering 
collisions we obtain the impurity transport coefficients 



(A17) 


Here lZz eS and Sz eff depend on the impurity content, for Z eS = 1.05 they are 77 i . 0 5 = 
-1.04549, <Si .05 = 1.82082 and for Z eS = 2.0 they are 77 2 . 0 = -1.14447, S 2 . 0 = 0.329224. 

The high-collisionality asymptotes for the W7-X case in Sec. 11 are plotted in Figs. 12 
and 13, and compared to SFINCS calculations at F* = 0 with both full linearized Fokker- 
Planck-Landau collisions and pitch-angle scattering (note that since SFINCS includes 
E x B-precession whereas the analytic high-collisionality calculations do not, a comparison 
at finite E * is not meaningful). We fold that the SFINCS results conform well with the 
analytic predictions. 
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♦ Z eil = 1.05, E* =0, SFINCS Fokker-Planck-Landau ■■■■ Z e{[ = 1.05, £*=0, FPL Pfirsch-Schluter asymptote 
•—• Z c([ = 1.05, E* =0, SFINCS Pitch-angle-scattering . Z c[[ = 1.05, E t =0, PAS Pfirsch-Schluter asymptote 



FIG. 12: Carbon (Z = 6) transport coefficients Lff (a), Lf\ (b), L\ 2 (c) and the temperature 
gradient coefficient Zf 2 — 3 + Lf^J /2 (d) as functions of normalized collision frequency u' z 

for a W7-X geometry at E* = 0 and with Z e g = 1.05. SFINCS computations for Fokker-Planck- 
Landau collisions ( —♦— ) and pitch-angle scattering (— •— ) are compared to the analytic high- 

collisionality limits for Fokker-Planck-Landau collisions ( - ) and pitch-angle scattering ( — ). 

Note the double-logarithmic scale in (b)-(d). 
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FIG. 13: Carbon (Z = 6) transport coefficients Lff (a), Ll\ (b), L\ 2 (c) and the temperature 
gradient coefficient L\ 2 ~ 3 + L\\ ^ /2 (d) as functions of normalized collision frequency u' z 

for a W7-X geometry at E* = 0 and with Z e g- = 2.0. SFINCS computations for Fokker-Planck- 
Landau collisions ( —♦— ) and pitch-angle scattering (— •— ) are compared to the analytic high- 

collisionality limits for Fokker-Planck-Landau collisions ( - ) and pitch-angle scattering ( — ). 

Note the double-logarithmic scale in (b)-(d). 
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